#define _BILUTEST_
program PRIDE
    use OMP_LIB
    use rng
    use cmfd2g,         only : initcmfd2g, solcmfd2g
    use lcmfd2g,        only : initlcmfd2g
    use cmfdmg,         only : initCMFDmg, initFDMMG, FDM, readFDM

#ifdef MPIF
	use MPIComm_Mod,    only : mpi_mc_comm
#endif

    use asmcal

    implicit none
#ifdef MPIF
	include "mpif.h"
#endif
    include "ntbytes.FI"
	include "xsec.FI"
    include "pspec.FI"
    real(8) ::k, err
    integer :: iter, ierr
    
#ifdef MPIF
	call mpi_init(ierr)
	call mpi_comm_size(mpi_mc_comm, nprocs, ierr)
	call mpi_comm_rank(mpi_mc_comm, myrank, ierr)
	call build_derive_type()
	allocate(nqmpi(0:nprocs-1))
	nqmpi=0
#ifdef MPIDEBUG
21  format ("NPROCS :", i2, 1x, "rank (", i2, ") is created.")
	print 21, nprocs, myrank
#endif
#endif
    call readinput()
    
#ifdef _THFDB
    call readxsec()
    call initTHFDB()
#endif    

!    call prepare_fdm()
!    call fdm(k,err,iter)            
    
#ifndef _NOCMFD
    call initcmfd2g
    call initlcmfd2g

!    call prepare_cmfd()
!    call cmfd(k, err, iter, .TRUE.)
#endif

#ifdef MPIF !=============================================
    if (myrank .eq. 0) then
#endif      !=============================================
    if (MC_ASMCAL) then
        tallytype = MD_ASMCAL
        call MC_ASM()
    endif
    call initcmfdmg 
    call initFDMMG()

    if (lmgfdmfirst .or. MC_ASMCAL) call FDM()
		
	if (lreadFDM) then
		print *, "Read initial FSD"
		call readFDM()
	endif

!    call initgeom
#ifdef MPIF !=============================================
    endif
#endif      !=============================================

!    call solcmfd2g
    
!    phif = 1.
!    call init_FMC
    print *, 1._8*nbytesf/1e6_8
    
!    call migrate_inasm
    
    tallytype = MD_MCCMFD
    call MC(1._8)
#ifdef MPIF
	call mpi_finalize(ierr)
#endif
end program

#ifdef MPIF
subroutine build_derive_type()
	use define
    use MPIComm_mod, only : mpi_mc_comm
	implicit none
	include "mpif.h"
	include "pspec.FI"
	integer :: i
	integer :: mesg_mpi_t
	integer :: block_lengths(3)
	integer :: ierr
	integer :: typelist(3)
	integer :: displacements(3)
	integer :: start_address, address
	type(particle) :: ntrn
	real(8) :: a, b, c
	type(particle), pointer, dimension(:) :: sn, rn
	integer, dimension(2) :: nqsend, rcvcnt, idisp
	allocate(sn(10), rn(10))
		
	block_lengths(1) = 2
	block_lengths(2) = 1
	block_lengths(3) = 2

	typelist(1) = MPI_DOUBLE_PRECISION
	typelist(2) = MPI_DOUBLE_PRECISION
	typelist(3) = MPI_INTEGER
	
	call MPI_ADDRESS(ntrn%xy,  start_address, ierr)
	displacements(1) = 0

	call MPI_ADDRESS(ntrn%w,  address, ierr)
	displacements(2) = address - start_address

	call MPI_ADDRESS(ntrn%lm,  address, ierr)
	displacements(3) = address - start_address

	call MPI_TYPE_STRUCT(3, block_lengths, displacements, &
			typelist, t_ntrn, ierr)
	call MPI_TYPE_COMMIT(t_ntrn, ierr)
!	print *, myrank, "Derived Type ID", t_ntrn
!	sn(1)%w=myrank*10 + 1.
!	sn(2)%w=myrank*10 + 2.
!	nqsend = (/0, 2 /)
!	rcvcnt = (/0, 2/)
!	idisp = (/0, 0/)
!	call mpi_allgatherv(sn(1), nqsend(myrank+1), t_ntrn, &
!						rn, rcvcnt, idisp, t_ntrn, &
!						mpi_mc_comm, ierr)
!	print *, myrank, rn(1)%w, rn(2)%w, rn(3)%w, rn(4)%w
!
end subroutine
#endif
